#Replication Code for IO Main Text (6.28.2023)

library(tidyverse)
library(skimr)
library(gplots)
library(broom)
library(dotwhisker)
library(lubridate)
library(lavaan)
library(jsonlite)
library(mediation)
library(psych)
library(ltm)
library(pwr)

## Quick fix for stargazer <= 5.2.3 is.na() issue with long model names in R >= 4.2
# Unload stargazer if loaded
detach("package:stargazer",unload=T)
# Delete it
remove.packages("stargazer")
# Download the source
download.file("https://cran.r-project.org/src/contrib/stargazer_5.2.3.tar.gz", destfile = "stargazer_5.2.3.tar.gz")
# Unpack
untar("stargazer_5.2.3.tar.gz")
# Read the sourcefile with .inside.bracket fun
stargazer_src <- readLines("stargazer/R/stargazer-internal.R")
# Move the length check 5 lines up so it precedes is.na(.)
stargazer_src[1990] <- stargazer_src[1995]
stargazer_src[1995] <- ""
# Save back
writeLines(stargazer_src, con="stargazer/R/stargazer-internal.R")
# Compile and install the patched package
install.packages("stargazer", repos = NULL, type="source")

library(stargazer)


# LOADING DATA
load(file = "final_data_IO.Rdata")

#SECTION 6 (RESULTS)

#6.1 EMOTIONS AND EXPOSURE TO TERRORISM NEWS

#Dominant emotional responses
final_data_IO  %>% filter(newsstory == "Attack") %>% group_by(angermore = anger > fear) %>% count(angermore)
final_data_IO  %>% filter(newsstory == "Attack") %>% group_by(fearmore = anger < fear) %>% count(fearmore)
final_data_IO  %>% filter(newsstory == "Attack") %>% group_by(equalemotions = anger == fear) %>% count(equalemotions)

#Anger higher across demographic categories
final_data_IO %>% group_by(ethnicity_cat) %>% summarize(mean(anger, na.rm = TRUE), mean(fear, na.rm = TRUE))
final_data_IO %>% group_by(male) %>% summarize(mean(anger, na.rm = TRUE), mean(fear, na.rm = TRUE))
final_data_IO %>% group_by(religion) %>% summarize(mean(anger, na.rm = TRUE), mean(fear, na.rm = TRUE))
final_data_IO %>% group_by(residence) %>% summarize(mean(anger, na.rm = TRUE), mean(fear, na.rm = TRUE))
final_data_IO %>% group_by(educ) %>% summarize(mean(anger, na.rm = TRUE), mean(fear, na.rm = TRUE))
final_data_IO %>% group_by(income) %>% summarize(mean(anger, na.rm = TRUE), mean(fear, na.rm = TRUE))

#Paired Sample T-Test comparing anger to fear across each condition
ttest_angerfear = subset(final_data_IO, nat_med == 1) %>% 
  dplyr::select(condition, anger, fear) %>%
  mutate(condition = case_when(
    condition == "True Control" ~ "No Attack",
    condition == "JC No Emotion" ~ "Jihadi-Civilian",
    condition == "JM No Emotion" ~ "Jihadi-Military",
    condition == "AC No Emotion" ~ "WN-Civilian",
    condition == "AM No Emotion" ~ "WN-Military"
  )) %>% 
  mutate(condition = fct_relevel(condition, "No Attack")) %>%
  pivot_longer(cols = -condition, names_to = "Emotion") %>% 
  mutate(Emotion = ifelse(Emotion == "fear", "Fear", "Anger")) %>% 
  group_by(condition) %>% 
  do(tidy(t.test(value~Emotion, paired = TRUE, data=.)))

#T-Test Table in Appendix
subset(final_data_IO, nat_med == 1) %>% 
  dplyr::select(condition, anger, fear) %>%
  pivot_longer(cols = c(anger, fear), names_to = "Emotion") %>% 
  group_by(condition) %>% 
  summarise("avgAnger" = mean(value[Emotion == "anger"], na.rm=TRUE), 
            "sdAnger" = sd(value[Emotion == "anger"], na.rm=TRUE), 
            "semAnger" = sd(value[Emotion == "anger"], na.rm=TRUE) / sqrt(n()),
            
            "avgFear" = mean(value[Emotion == "fear"], na.rm=TRUE), 
            "sdFear" = sd(value[Emotion == "fear"], na.rm=TRUE), 
            "semFear" = sd(value[Emotion == "fear"], na.rm=TRUE) / sqrt(n()),
            
            t.test_stat = t.test(paired = TRUE, value[Emotion == "fear"], value[Emotion == "anger"])$statistic,
            t.test_pval = t.test(paired = TRUE, value[Emotion == "fear"], value[Emotion == "anger"])$p.value,
            t.test_stderr = t.test(paired = TRUE, value[Emotion == "fear"], value[Emotion == "anger"])$stderr
  ) %>% 
  knitr::kable(format = "latex")

#Figure 3a: Relative levels of anger and fear in response to terrorism news

#pdf(width = 8, height = 5, file = "barplot_emotions.pdf")

subset(final_data_IO, nat_med == "1") %>%
  dplyr::select(condition, anger, fear) %>%
  mutate(condition = case_when(
    condition == "True Control" ~ "No Attack",
    condition == "JC No Emotion" ~ "Jihadi-Civilian",
    condition == "JM No Emotion" ~ "Jihadi-Military",
    condition == "AC No Emotion" ~ "WN-Civilian",
    condition == "AM No Emotion" ~ "WN-Military"
  )) %>%
  mutate(condition = fct_relevel(condition, "No Attack")) %>%
  pivot_longer(cols = -condition, names_to = "Emotion") %>%
  mutate(Emotion = ifelse(Emotion == "fear", "Fear", "Anger")) %>%
  ggplot(aes(x = condition, y = value, fill = Emotion)) +
  stat_summary(geom = "bar", position = "dodge") +
  stat_summary(geom = "errorbar", position = position_dodge(0.95), width = 0.3) +
  scale_fill_manual(name = "Emotion", values = RColorBrewer::brewer.pal(11,"Spectral")[c(2, 11)]) +
  coord_cartesian(ylim = c(1,5)) +
  labs(x = "Attack Types", y = "Self Reported Emotions") +
  cowplot::theme_cowplot()
# scale_fill_grey()

#dev.off()


#Figure 3b: Relative levels of anger and fear in response to terrorism news
reg_ATEanger <- lm(anger ~ condition, data=subset(final_data_IO, nat_med=="1"))
summary(reg_ATEanger)

reg_ATEfear <- lm(fear ~ condition, data=subset(final_data_IO, nat_med=="1"))
summary(reg_ATEfear)

ATEemotions_model <-
  bind_rows("Anger" = tidy(reg_ATEanger),
            "Fear" = tidy(reg_ATEfear),
            .id = "model") %>% 
  filter(term != "(Intercept)") %>% 
  mutate(model = (fct_inorder(model)),
         term = case_when(term == "conditionJC No Emotion" ~ "Jihadi-Civilian",
                          term == "conditionJM No Emotion" ~ "Jihadi-Military",
                          term == "conditionAC No Emotion" ~ "WN-Civilian",
                          term == "conditionAM No Emotion" ~ "WN-Military"
         )) %>% 
  mutate(conf.low90 = estimate - stats::qnorm(1 - (1 - .90) /
                                                            2) * std.error,
        conf.high90 = estimate + stats::qnorm(1 - (1 - 0.90) /
                                                             2) * std.error)

#pdf(width = 7, height = 5, file = "ATEemotions_model.pdf")

dwplot(ATEemotions_model,ci = 0.95, 
       dot_args = list(size = 4, aes(shape = model, color = model)),
       whisker_args = list(size = 1, aes(color = model))) +
  ggstance::geom_linerangeh(aes(
                           y = term,
                           xmin = conf.low90,
                           xmax = conf.high90,
                           group = interaction(model, term),
                           color = model),
                           size = 2,
                           position = ggstance::position_dodgev(height = 0.4),
                           ) + 
  geom_vline(xintercept = 0, linetype = 2) + 
  scale_color_manual(name = "Emotion", values = RColorBrewer::brewer.pal(11,"Spectral")[c(11, 2)]) +
  # scale_shape_manual(name = "Emotion", values = c(16, 17)) +
  guides(color = guide_legend(reverse=T)) +
  labs(y = "Attack Type", x = "Effect of Exposure to Attack on Self-Reported Emotions") +
  cowplot::theme_cowplot() +
  guides(shape = guide_legend("Emotion"), colour = guide_legend("Emotion"))
  # scale_color_grey()

#dev.off()

#Respondents who live geographically close also exhibit more anger than fear (Footnote 67) 
final_data_IO  %>% count(WyomingNeighbor)

final_data_IO %>% 
  group_by(newsstory, WyomingNeighbor) %>% 
  summarize(mean(anger, na.rm = TRUE), mean(fear, na.rm = TRUE))



#SECTION 6.2 Effects of Exposure to Terrorism News on Political Attitudes

#ATE of News Exposure on Attitudes
reg_ATEretaliate <- lm(retaliate ~ newsstory, data=filter(final_data_IO, nat_med =="1"))
summary(reg_ATEretaliate)

reg_ATEconciliate <- lm(conciliate ~ newsstory, data=subset(final_data_IO, nat_med=="1"))
summary(reg_ATEconciliate)

reg_ATEpunitive <- lm(punish ~ newsstory, data=filter(final_data_IO, nat_med =="1"))
summary(reg_ATEpunitive)


#Table 1: Support for retaliation by terrorism type and partisanship

#model 1
summary(reg_ATEretaliate)

#model 2
reg_ATEretaliate_victimid <- lm(retaliate ~ victim_id, data=filter(final_data_IO, nat_med =="1")%>% 
                                  mutate(victim_id = fct_relevel(victim_id, "Civilian")))
summary(reg_ATEretaliate_victimid)

#model 3
reg_ATEretaliate_attackerid <- lm(retaliate ~ attacker_id, data=filter(final_data_IO, nat_med =="1")%>% 
                                    mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")))
summary(reg_ATEretaliate_attackerid)

#model 4
reg_INTretaliate_attackerid <- lm(retaliate ~ attacker_id*pid_group_merged, data=subset(final_data_IO, nat_med=="1")%>% 
                                    mutate(pid_group_merged = fct_relevel(pid_group_merged, "Democrat"))%>% 
                                    mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")))
summary(reg_INTretaliate_attackerid)

#model 5
reg_INTretaliate_victimid <- lm(retaliate ~ victim_id*pid_group_merged, data=subset(final_data_IO, nat_med=="1")%>% 
                                  mutate(pid_group_merged = fct_relevel(pid_group_merged, "Democrat"))%>% 
                                  mutate(victim_id = fct_relevel(victim_id, "Civilian")))
summary(reg_INTretaliate_victimid)

stargazer(reg_ATEretaliate, reg_ATEretaliate_victimid, reg_ATEretaliate_attackerid, reg_INTretaliate_attackerid, reg_INTretaliate_victimid, align=TRUE)


#Figure 4: Interaction of partisanship and perpetrator ID on responses to terorrism news (comparison to Jihadist)

#pdf(width = 8, height = 5, file = "INTangerperpID.pdf")

final_data_IO %>% 
  filter(nat_med == "1") %>% 
  filter(!is.na(pid_group_merged)) %>% 
  mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")) %>%   
  group_by(pid_group_merged) %>% 
  do(broom::tidy(lm(anger ~ attacker_id, data = .))) %>% 
  rename(model = pid_group_merged) %>% 
  relabel_predictors(
    c(attacker_idAryan = "WN Attacker")
  ) %>% 
  filter(term != "attacker_idNone") %>% 
  # filter(model != "Independent") %>% 
  bind_rows("Self-Reported Anger" = .,
            "Punitive Motives" = 
              final_data_IO %>% 
              filter(!is.na(pid_group_merged)) %>% 
              filter(nat_med == "1") %>% 
              mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")) %>%   
              group_by(pid_group_merged) %>% 
              do(broom::tidy(lm(punish ~ attacker_id, data = .))) %>% 
              rename(model = pid_group_merged) %>% 
              relabel_predictors(
                c(attacker_idAryan = "WN Attacker")
              ) %>% 
              filter(term != "attacker_idNone"), # %>% 
              # filter(model != "Independent"),
            "Support for Retaliation" = 
              final_data_IO %>% 
              filter(!is.na(pid_group_merged)) %>% 
              filter(nat_med == "1") %>% 
              mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")) %>%   
              group_by(pid_group_merged) %>% 
              do(broom::tidy(lm(retaliate ~ attacker_id, data = .))) %>% 
              rename(model = pid_group_merged) %>% 
              relabel_predictors(
                c(attacker_idAryan = "WN Attacker")
              ) %>% 
              filter(term != "attacker_idNone"), # %>% 
            # filter(model != "Independent"),
            .id = "DV") %>% 
  mutate(DV = fct_relevel(DV, "Self-Reported Anger")) %>%
  mutate(conf.low90 = estimate - stats::qnorm(1 - (1 - .90) /
                                                2) * std.error,
         conf.high90 = estimate + stats::qnorm(1 - (1 - 0.90) /
                                                 2) * std.error) %>% 
  dwplot(dot_args = list(size = 4, aes(shape = model, color = model)),
         whisker_args = list(size = 1, aes(color = model))) +
  ggstance::geom_linerangeh(aes(
    y = term,
    xmin = conf.low90,
    xmax = conf.high90,
    group = interaction(model, term),
    color = model),
    size = 2,
    position = ggstance::position_dodgev(height = 0.4),
  ) + 
  geom_vline(xintercept = 0, linetype = 2) +
  scale_color_manual(name = "Party ID", values = RColorBrewer::brewer.pal(3,"Set1")[c(1, 3, 2)]) +
  guides(color = guide_legend(reverse=T)) +
  labs(x = "Effect of White Nationalist Attack Compared to Jihadist Attack") +
  cowplot::theme_cowplot() +
  cowplot::panel_border() +
  guides(shape = guide_legend("Party ID"), colour = guide_legend("Party ID")) +
  facet_grid(~ DV)

#dev.off()


#SECTION 6.3 The Causal Impact of Anger

#Figure 5: Effect of Distinct Emotions after Exposure to Terrorism on Policy Preferences
reg_ANIEretaliate_noprime <- lm(retaliate ~ emotion, data=subset(final_data_IO %>% filter(newsstory != "None")))
summary(reg_ANIEretaliate_noprime)

reg_ANIEconciliate_noprime <- lm(conciliate ~ emotion, data = subset(final_data_IO %>% filter(newsstory != "None")))
summary(reg_ANIEconciliate_noprime)


#pdf(width = 7, height = 5, file = "ANIE_policies.pdf")

ANIEpolicy_model <-
  bind_rows("Retaliate" = tidy(reg_ANIEretaliate_noprime),
            "Conciliate" = tidy(reg_ANIEconciliate_noprime),
            .id = "model") %>% 
  filter(term != "(Intercept)", term != "emotionNone", term != "emotionReassure") %>% 
  mutate(model = (fct_inorder(model)),
         term = str_remove(term, "emotion")) %>% 
  mutate(conf.low90 = estimate - stats::qnorm(1 - (1 - .90) /
                                                2) * std.error,
         conf.high90 = estimate + stats::qnorm(1 - (1 - 0.90) /
                                                 2) * std.error)

dwplot(ANIEpolicy_model %>% mutate(model = fct_relevel(model, "Conciliate")), 
       dot_args = list(size = 4, aes(shape = model, color = model)),
       whisker_args = list(size = 1, aes(color = model))) +
  ggstance::geom_linerangeh(aes(
    y = term,
    xmin = conf.low90,
    xmax = conf.high90,
    group = interaction(model, term),
    color = model),
    size = 2,
    position = ggstance::position_dodgev(height = 0.4),
  ) + 
  geom_vline(xintercept = 0, linetype = 2) + 
  scale_color_manual(name = "Policy", values = RColorBrewer::brewer.pal(11,"RdYlGn")[c(8, 11)]) +
  guides(color = guide_legend(reverse=T)) +
  labs(y = "Emotion Prime", x = "Effect of Emotion Prime on Political Attitudes") +
  cowplot::theme_cowplot() +
  # scale_color_grey(name = "Policy") +
  guides(shape = guide_legend("Policy", reverse = T), colour = guide_legend("Policy", reverse = T))

#dev.off()



#Figure 6: Effect of Distinct Emotions after Exposure to Terrorism News on Motives for Retaliation
reg_ANIErehab <- lm(rehab ~ emotion, data = subset(final_data_IO %>% filter(newsstory != "None")))
summary(reg_ANIErehab)

reg_ANIEincapac <- lm(incapac ~ emotion, data = subset(final_data_IO %>% filter(newsstory != "None")))
summary(reg_ANIEincapac)

reg_ANIEdeter <- lm(deter ~ emotion, data = subset(final_data_IO %>% filter(newsstory != "None")))
summary(reg_ANIEdeter)

reg_ANIEpunish <- lm(punish ~ emotion, data = subset(final_data_IO %>% filter(newsstory != "None")))
summary(reg_ANIEpunish)

ANIE_models <-
  bind_rows("Rehabilitate" = tidy(reg_ANIErehab),
            "Incapacitate" = tidy(reg_ANIEincapac),
            "Deter" = tidy(reg_ANIEdeter),
            "Punish" = tidy(reg_ANIEpunish),
            .id = "model") %>% 
  filter(term != "(Intercept)", term != "emotionNone", term != "emotionReassure") %>% 
  mutate(model = (fct_inorder(model)),
         term = str_remove(term, "emotion")) %>% 
  mutate(conf.low90 = estimate - stats::qnorm(1 - (1 - .90) /
                                                2) * std.error,
         conf.high90 = estimate + stats::qnorm(1 - (1 - 0.90) /
                                                 2) * std.error)

#pdf(width = 7, height = 5, file = "ANIE_motives.pdf")

dwplot(ANIE_models %>% mutate(model = (fct_relevel(model, "Punish"))), 
       dot_args = list(size = 4, aes(shape = model, color = model)),
       model_order = c("Punish", "Deter", "Incapacitate", "Rehabilitate"),
       whisker_args = list(size = 1, aes(color = model))) +
  ggstance::geom_linerangeh(aes(
    y = term,
    xmin = conf.low90,
    xmax = conf.high90,
    group = interaction(model, term),
    color = model),
    size = 2,
    position = ggstance::position_dodgev(height = 0.4),
  ) + 
  geom_vline(xintercept = 0, linetype = 2) + 
  scale_color_manual(name = "Motive", values = RColorBrewer::brewer.pal(11,"Spectral")[c(5, 2, 3, 1)]) +
  guides(color = guide_legend(reverse=T)) +
  labs(y = "Emotion Prime", x = "Effect of Emotion Prime on Motives for Retaliation") +
  cowplot::theme_cowplot() +
  # scale_color_grey(name = "Motive")+
  guides(shape = guide_legend("Motive", reverse = T), colour = guide_legend("Motive", reverse = T))


#dev.off()

#Anger and fear coefficients not significantly different from each other (Footnote 72)
reg_ANIEpunish_fear <- lm(punish ~ emotion, data = subset(final_data_IO %>% filter(newsstory != "None") %>%
                                                            mutate(emotion = fct_relevel(emotion, "Fear")))
)
summary(reg_ANIEpunish_fear)


#Figure 7: Effect of emotions primes on the interaction of partisanship and perpetrator ID on attitudes

#Interaction effect of PIDxPerpetrator with NO PRIME
reg_ACDE_retaliate <- lm(retaliate ~ pid_group_merged*attacker_id, data = subset(final_data_IO, emotion == "None")%>% 
                           mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")))
summary(reg_ACDE_retaliate)

#Interaction effect of PIDxPerpetrator with ANGER PRIME
reg_ACDE_retaliate1 <- lm(retaliate ~ pid_group_merged*attacker_id+RWA_normal+NFC_normal+trait_anger_normal, data = subset(final_data_IO, emotion == "Anger")%>% 
                            mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")))
summary(reg_ACDE_retaliate1)

#Interaction effect of PIDxPerpetrator with NO PRIME
reg_ACDE_punish <- lm(punish ~ pid_group_merged*attacker_id, data = subset(final_data_IO, emotion == "None")%>% 
                        mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")))
summary(reg_ACDE_punish)

#Interaction effect of PIDxPerpetrator with ANGER PRIME
reg_ACDE_punish1 <- lm(punish ~ pid_group_merged*attacker_id, data = subset(final_data_IO, emotion == "Anger")%>% 
                         mutate(attacker_id = fct_relevel(attacker_id, "Jihadist")))
summary(reg_ACDE_punish1)

#pdf(width = 8, height = 5, file = "ACDEattitudes_PID.pdf")

final_data_IO %>% 
  filter(attacker_id == "Jihadist" | attacker_id == "Aryan") %>% 
  filter(emotion == "Anger" | emotion == "None") %>% 
  filter(pid_group_merged == "Democrat" | pid_group_merged == "Republican") %>% 
  dplyr::select(attacker_id, emotion, punish,  retaliate, pid_group_merged) %>% 
  pivot_longer(cols = c(punish, retaliate), names_to = "DVs") %>%
  # mutate(emotion = fct_relevel(emotion, "None", "Anger", "Fear")) %>% 
  mutate(DVs = case_when(DVs == "punish" ~ "Punitive motives",
                         DVs == "retaliate" ~ "Support for retaliation")) %>% 
  ggplot(aes(x = attacker_id, y = value, group = pid_group_merged, color = pid_group_merged, shape = pid_group_merged)) + 
  stat_summary(fun.data = mean_cl_normal) +
  stat_summary(geom = "line") +
  facet_grid(DVs ~ emotion) +
  scale_color_manual(name = "Party ID", values = RColorBrewer::brewer.pal(2,"Set1")[c(2, 1)]) +
  scale_x_discrete(labels = c("Jihadist", "White Nationalist")) +
  labs(x = "Perpetrator", y = "Attitude") +
  coord_cartesian(ylim = c(5, 7)) +
  cowplot::theme_cowplot() +
  cowplot::panel_border() +
  guides(shape = guide_legend("Party ID"), colour = guide_legend("Party ID"))

#dev.off()


#No interaction of AEMT response by partisanship and perpetrator ID (triple interaction) (Footnote 75)

reg_ANIEangerxPID2 <- lm(anger ~ emotion*pid_group_merged*attacker_id, data=filter(final_data_IO, newsstory!="None"))
summary(reg_ANIEangerxPID2)

reg_ANIEfearxPID2 <- lm(fear ~ emotion*pid_group_merged*attacker_id, data=filter(final_data_IO, newsstory!="None"))
summary(reg_ANIEfearxPID2)



#SECTION 7: ANGER AND FEAR AFTER REAL WORLD ATTACKS

load(file = "twitter_data_IO.Rdata")

tweet_data %>% 
  dplyr::select(Anxiety, Anger, Anxiety.Risk, Anger.Morality,attack, attackLabel, attackDate) %>% 
  pivot_longer(cols = c(Anxiety, Anger, Anxiety.Risk, Anger.Morality), names_to = "Content") %>% 
  group_by(attack, Content) %>% 
  summarise("avg" = mean(value), sd(value), "sem" = sd(value) / sqrt(n()))


#Figure 8: Sentiment Analysis of Tweets Surrounding Real World Terror Attacks
#pdf(width = 10, height = 5, file = "tweets_bar_ci2.pdf")

tweet_data %>% 
  dplyr::select(Anxiety, Anger, attack, attackLabel, attackDate) %>% 
  pivot_longer(cols = c(Anxiety, Anger), names_to = "Emotion") %>% 
  mutate(Date = lubridate::ymd(attackDate)) %>% 
  mutate(attackLabel = fct_reorder(attackLabel, .x = attackDate, .fun = max)) %>% 
  ggplot(aes(x = attackLabel, y = value, fill = Emotion)) +
  stat_summary(geom = "bar", position = "dodge") +
  stat_summary(width = 0.2, 
               geom = "errorbar", 
               position = position_dodge(0.95), 
               fun.data = mean_cl_normal) +
  labs(y = "Avg. Tweet Sentiment %", x = "Attack") +
  scale_fill_manual(name = "Emotion", 
                    labels = c("Anger", "Anxiety"), 
                    values = RColorBrewer::brewer.pal(11,"Spectral")[c(2, 11)]) +
  cowplot::theme_cowplot() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

#dev.off() 


#T-Test Analysis
tweet_data %>% 
  dplyr::select(Anxiety, Anger, attack, attackDate) %>% 
  pivot_longer(cols = c(Anxiety, Anger), names_to = "Emotion") %>% 
  group_by(attack) %>% 
  mutate(Date = lubridate::ymd(attackDate),
         Attack = attack) %>% 
  summarise("Mean Anger" = mean(value[Emotion == "Anger"]), 
            "SD Anger" = sd(value[Emotion == "Anger"]), 
            #"semAnger" = sd(value[Emotion == "Anger"]) / sqrt(n()),
            
            "Mean Fear" = mean(value[Emotion == "Anxiety"]), 
            "SD Fear" = sd(value[Emotion == "Anxiety"]), 
            #"semFear" = sd(value[Emotion == "Anxiety"]) / sqrt(n()),
            
            t.test_stat = t.test(paired = TRUE, value[Emotion == "Anxiety"], value[Emotion == "Anger"])$statistic,
            t.test_df = t.test(paired = TRUE, value[Emotion == "Anxiety"], value[Emotion == "Anger"])$parameter,
            t.test_pval = t.test(paired = TRUE, value[Emotion == "Anxiety"], value[Emotion == "Anger"])$p.value,
            #t.test_stderr = t.test(paired = TRUE, value[Emotion == "Anxiety"], value[Emotion == "Anger"])$stderr
            Date = max(attackDate)
  ) %>% 
  arrange(Date) %>%
  knitr::kable(format = "latex", digits = 2)
